Цепочка Вольтерра, численный счёт

В.Э. Адлер. Классические интегрируемые системы, весенний семестр 2020

Несколько простых численных экспериментов с цепочкой Вольтерра

(Subscript[u, n,t])=(Subscript[u, n]) ((Subscript[u, n+1])-(Subscript[u, n-1])).

Решение солитонного типа

Активируйте ячейки нажатием [Shift-Enter].

(* задаём период и время счёта *)
M1=100;
t1=1000;

(* уравнения, с учётом периодического замыкания с периодом M *)
eqs=Table[u[n]'[t]==u[n][t](u[Mod[n+1,M1]][t]-u[Mod[n-1,M1]][t]),{n,0,M1-1}];

(* задаём начальные условия. Попробуйте что-нибудь поменять *)
ic=Table[{u[n][0]==1+0.5Exp[-0.01(((n-M1/2))^2)]},{n,0,M1-1}];

(* строим решение при помощи NDSolve *)
sol1=NDSolve[Join[eqs,ic],Table[u[n],{n,0,M1-1}],{t,0,t1},MaxSteps->∞][[1]];

(* массив точек *)
u1[t_]:=Table[{n,u[n][t]},{n,0,M1-1}]/.sol1

(* функция для вывода. Для наглядности, массив выводится два раза: как ломаная линия 
и как набор жирных точек в ее вершинах *)
gr1[t_]:=ListPlot[{u1[t],u1[t]},
PlotRange->{{-0.1,M1+1},{0.95,2.05}},
PlotStyle->{{Thin,Black},{Blue,PointSize[0.01]}},
Joined->{True,False},
AxesLabel->{" \!\(\*SubscriptBox[\(n\), \(\\ \)]\)","(Subscript[u, n])"},
BaseStyle->{FontSize->14,FontFamily->"Times New Roman"},
ImageSize->400,
Epilog->{Text["\!\(\*SubscriptBox[\(t\), \(\\ \)]\)= "<>ToString[t],{0.8M1,1.95},{-1,0}]}
]

(* анимация *)
Manipulate[gr1[t],{{t,0},0,t1}]

Имитация эксперимента Забуски-Краскала

Если брать начальное условие с более широким горбиком, то мы придем примерно к той картине, что была в эксперименте Забуски-Краскала. Напомним, что там решалось чисто разностное уравнение с начальными условиями в виде синусоидальной волны с большим периодом (порядка 100 узлов решетки) и не очень большой амплитудой. Такое начальное условие распадалось на несколько солитонов плюс “рябь” (число солитонов зависит от начальной амплитуды), которые двигались с разной скоростью и, благодаря периодичности, всё время нагоняли друг друга. Качественно, в цепочке Вольтерра наблюдается точно такая же картина.

M2=150;
t2=3000;

eqs=Table[u[n]'[t]==u[n][t](u[Mod[n+1,M2]][t]-u[Mod[n-1,M2]][t]),{n,0,M2-1}];
ic=Table[{u[n][0]==1+0.07Cos[((2π)/M2) n]},{n,0,M2-1}];
sol2=NDSolve[Join[eqs,ic],Table[u[n],{n,0,M2-1}],{t,0,t2},MaxSteps->∞][[1]];

u2[t_]:=Table[{n,u[n][t]},{n,0,M2-1}]/.sol2

gr2[t_]:=ListPlot[{u2[t],u2[t]},
PlotRange->{{-0.1,M2+1},{0.9,1.3}},
PlotStyle->{{Thin,Black},{Blue,PointSize[0.007]}},
Joined->{True,False},
AxesLabel->{" \!\(\*SubscriptBox[\(n\), \(\\ \)]\)","(Subscript[u, n])"},
BaseStyle->{FontSize->14,FontFamily->"Times New Roman"},
ImageSize->400,
Epilog->{Text["\!\(\*SubscriptBox[\(t\), \(\\ \)]\)= "<>ToString[t],{0.8M2,1.25},{-1,0}]}
]

Manipulate[gr2[t],{{t,0},0,t2}]

Игра с чётностью

В цепочке Вольтерра есть стационарное решение (Subscript[u, 2n])=a,  (Subscript[u, 2n+1])=b с разными постоянными a и b. Если немного возмутить его, взяв начальное условие в виде двух плавных кривых, отдельно для чётных и нечётных узлов, то эти узлы так и останутся разделёнными. Конечно, в уравнении КдФ такое невозможно. Для решений такого типа тоже можно написать непрерывный предел, но уже к системе на две полевые переменные.

В этом примере поведение решения кажется более сложным и счёт не такой устойчивый; если сильно увеличить амплитуду возмущений, то в решении могут возникнуть особенности.

M3=200;
t3=10000;

eqs=Table[u[n]'[t]==u[n][t](u[Mod[n+1,M3]][t]-u[Mod[n-1,M3]][t]),{n,0,M3-1}];

ic=Table[{u[n][0]==
If[OddQ[n],1+0.1Exp[-0.005(((n-M3/3))^2)],   (* формула для нечётных узлов *) 
1.2-0.05Exp[-0.005(((n-2M3/3))^2)]]},            (* формула для чётных узлов *)
{n,0,M3-1}];

sol3=NDSolve[Join[eqs,ic],Table[u[n],{n,0,M3-1}],{t,0,t3},MaxSteps->∞][[1]];

u3[t_]:=Table[{n,u[n][t]},{n,0,M3-1}]/.sol3
u3o[t_]:=Table[{n,u[n][t]},{n,1,M3-1,2}]/.sol3
u3e[t_]:=Table[{n,u[n][t]},{n,0,M3-1,2}]/.sol3

gr3[t_]:=ListPlot[{u3[t],u3o[t],u3e[t],u3o[t],u3e[t]},
PlotRange->{{-0.1,M3+1},{0.85,1.55}},
PlotStyle->{{Thin,Black},{Thin,Black},{Thin,Black},
{Blue,PointSize[0.01]},{Red,PointSize[0.01]}},
Joined->{True,True,True,False,False},
AxesLabel->{" \!\(\*SubscriptBox[\(n\), \(\\ \)]\)","(Subscript[u, n])"},
BaseStyle->{FontSize->14,FontFamily->"Times New Roman"},
ImageSize->400,
Epilog->{Text["\!\(\*SubscriptBox[\(t\), \(\\ \)]\)= "<>ToString[t],{0.8M3,1.50},{-1,0}]}
]

Manipulate[gr3[t],{{t,0},0,t3}]

Упражнение

Это просто легкое упражнение на программирование, только для тех, кто осваивает язык.

Определите функцию (Subscript[u, n])(t) по явным формулам

(Subscript[u, n])(t)=(((Subscript[w, n+1])(t)(Subscript[w, n-2])(t))/((Subscript[w, n])(t)(Subscript[w, n-1])(t))),   (Subscript[w, n])(t)=1+c (a^n)exp((a-(a^(-1)))x)

где a и c произвольные постоянные, и проверьте, что она удовлетворяет цепочке Вольтерра. Постройте график, лучше в манипуляторе, чтобы можно было менять t и параметры.

То же самое для

(Subscript[u, n])=-((b(n+d)(e^(2b t)))/(a + (e^(2b t)))) при чётных n и (Subscript[u, n])=((b (a n + c))/(a + (e^(2b t)))) при нечётных n.

Это решение растет по n (поэтому при построении графика нужно выбрать масштаб покрупнее), то есть, оно не относится к солитонному типу. Но, тем не менее, это некоторое явное решение.